
// STL
#include <vector>     // std::vector<>

// jScience
#include "jstats.hpp" // mean() << jmm: temp >>

// stats++
#include "statsxx/statistics.hpp" // r


//========================================================================
//========================================================================
//
// NAME: double statistics::r2(
//                             const std::vector<double> &x,
//                             const std::vector<double> &y
//                             )
//
// DESC: Computes the square of the correlation coefficient (coefficient of determination).
//
// NOTES:
//
//     - See:
//
//     http://mathworld.wolfram.com/CorrelationCoefficient.html
//
//     https://en.wikipedia.org/wiki/Coefficient_of_determination#As_squared_correlation_coefficient
//
//     - r^2 is the proportion of ss_yy (see below) that is accounted for by the regression.
//
//========================================================================
//========================================================================
inline double statistics::r2(
                             const std::vector<double> &x,
                             const std::vector<double> &y
                             )
{
    // NOTE: I am not sure which of the following is faster, but the latter (calling r()) uses less code
    
/*
    auto n = x.size();
    
    double x_mean = mean(x);
    double y_mean = mean(y);
    
    double ss_xx = 0.0;
    double ss_yy = 0.0;
    double ss_xy = 0.0;
    
    for(auto i = 0; i < n; ++i)
    {
        ss_xx += (x[i]*x[i]);
        ss_yy += (y[i]*y[i]);
        ss_xy += (x[i]*y[i]);
    }
    
    ss_xx -= n*x_mean*x_mean;
    ss_yy -= n*y_mean*y_mean;
    ss_xy -= n*x_mean*y_mean;
    
    return (ss_xy*ss_xy/(ss_xx*ss_yy));
*/
    
    double r = statistics::r(x,y);
    
    return (r*r);
}
